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ABSTRACT 

A new  one-dimensional  bulk  model  of  the  mixed  layer  of  the  upper  ocean  is  presented.  An  entrainment 
hypothesis  dependent  upon  the  relative  distribution  of  turbulent  energy  between  horizontal  and  vertical 
components  is  offered  as  a plausible  mechanism  for  governing  both  entrainment  and  layer  retreat. 

This  model  has  two  properties  not  previously  demonstrated : 

(i)  The  fraction  of  wind-generated  turbulent  kinetic  energy  partitioned  to  potential  energy  increase 
by  means  of  mixed  layer  deepening  is  dependent  upon  layer  stability,  H*  = h/L,  as  measured  by  the  ratio 
of  mixed  layer  depth  h to  Obukhov  length  L.  This  results  in  a modulation  of  the  mean  entrainment  rate 
by  the  diurnal  heating  and  cooling  cycle. 

(ii)  Viscous  dissipation  is  enhanced  for  increased  values  of  Ro"‘  = A//w«,  where  / is  the  Coriolis  parame- 
ter and  u,  the  friction  velocity  for  the  water.  This  enables  a cyclical  steady  state  to  occur  over  an  annual 
period  by  limiting  maximum  layer  depth. 

A nondimensional  framework  used  to  present  the  general  solution  also  suggests  a basis  for  model  com- 
parison and  data  analysis. 
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1.  Introduction 

The  ocean  mixed  layer  treated  here  is  that  fully 
turbulent  region  of  the  upper  ocean  that  is  bounded 
above  by  the  sea-air  interface.  The  wind  and  intermit- 
tent upward  surface  buoyancy  flux  through  the  surface 
are  the  sources  of  mechanical  energy  for  the  generation 
of  this  turbulence.  Typically,  the  mixed  layer  is  bounded 
below  by  a dynamically  stable  watermass.  The  vertical 
fluxes  of  heat,  salt  and  momentum  in  the  turbulent 
boundary  layer  or  mixed  layer  are  essentially  decoupled 
from  those  of  the  underlying  stable  water  column  be- 
cause the  energy  for  this  mixing  comes  from  above. 
Minimal  vertical  fluxes  below  the  mixed  layer,  together 
with  high  turbulence  intensity  within  the  layer,  result 
in  an  approximate  vertical  uniformity  in  mean  velocity 
and  density.  This  ostensible  homogeneity  is  the  root  of 
the  term  “slab,”  which  is  often  used  to  describe  the 
layer  as  depicted  in  Fig.  1.  There  is  an  appealing 
practical  aspect  to  the  judicious  use  of  the  assumption 
of  vertical  homogeneity  in  a bulk  model  such  as  this 
because  the  problem  of  solving  for  the  interior  fluxes  of 
buoyancy  and  momentum  is  reduced  to  the  need  to 
know  only  the  surface  and  entrainment  fluxes.  However, 
only  small  vertical  gradients  in  these  mean  variables 
may  be  associated  with  large  turbulent  fluxes.  Therefore 


1 The  author  was  supported  during  this  research  at  the  Uni- 
versity of  Washington,  Seattle,  by  NOAA's  Environmental 
Research  Laboratories  and  by  NOAA's  GATE  Office. 


the  slab  assumption  should  not  be  as  readily  applied  to 
the  turbulent  kinetic  energy  budget. 

Earlier  works  of  concern  here  are  those  dealing 
explicitly  with  equations  for  the  production,  alteration 
and  destruction  of  turbulent  kinetic  energy  within  the 
mixed  layer.  Kraus  and  Turner  (1%7)  were  the  first 
to  heed  the  turbulent  kinetic  energy  budget  in  the 
prototype  one-dimensional  mixed  layer  model,  utilizing 
the  approximately  decoupled  state  of  the  equations  for 
the  thermal  and  mechanical  energies.  By  neglecting 
the  frictional  generation  of  heat,  the  vertically  inte- 
grated heat  equation  provides  a relationship  for  the 
conservation  of  potential  energy.  However,  viscous 
dissipation  cannot  be  neglected  in  the  turbulent  kinetic 
energy  budget.  This  last  aspect  has  been  recognized 
only  more  recently.  Dissipation  has  been  assumed  to  be 
a fixed  fraction  of  wind  stress  production  in  the  models 
of  Geisler  and  Kraus  (1969),  Miropol’skiy  (1970), 
Denman  (1973)  and  Niiler  (1975),  all  variations  of  the 
Kraus-Turner  model.  The  latest  parameterizations  of 
dissipation,  those  by  Elsberry  et  al.  (1976),  Resnyanskiy 
(1975)  and  Kim  (1976),  have  been  conceived  with  the 
recognition  of  a need  to  augment  dissipation  in  certain 
instances. 

These  earlier  theories  that  are  based  on  the  turbulent 
energy  equation  demonstrate  the  importance  of  the 
Obukhov  length  scale  I.  [as  first  applied  to  the  ocean 
by  Kitaigorodsky  (I960)].  However,  they  have  not 
explained  the  significance  of  another  length  scale, 
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Fig.  1.  Idealized  density  and  mean  velocity  profiles  of  the  ocean  mixed  layer. 
Region  I is  the  fully-turbulent  (Rf<Rf„)  mixed  layer  of  depth  h.  Region  II  is  the 
slightly  stable  (Rf=R3~l),  intermittently-turbulent  entrainment  zone  of  thickness 
i.  Region  III  is  the  stable  (Ri>Ri„)  underlying  watermass  having  negligible 
vertical  fluxes  in  comparison  to  those  of  region  I.  Rf  and  Ri  are  the  flux  and  gradient 
Richardson  numbers,  respectively. 


«,//,  as  proposed  by  Rossby  and  Montgomery  (1935). 
Specifically,  the  general  applicability  of  these  models  is 
limited  by  the  following  three  problems: 

1)  Gill  and  Turner  have  demonstrated  the  inabil- 
ity of  the  prototype  model  to  predict  cyclical  steady 
state  on  an  annual  basis.  If  the  viscous  dissipation  of 
turbulent  energy  is  parameterized  as  a fixed  fraction  of 
production,  the  fraction  of  wind-generated  energy  going 
to  entrainment  mixing  is  constrained  to  be  constant  by 
the  integral  relationship,  regardless  of  the  layer  depth. 
Testing  such  a model  Camp  (1976)  shows  that  the 
predicted  oeepening  during  storms  is  much  too  rapid 
and  unchecked.  An  enhancement  of  dissipation  is  one 
possible  answer,  but  a physical  explanation  is  required. 
In  a departure  from  the  Kraus-Turner  approach, 
Pollard  et  al.  (1973)  used  the  total  kinetic  energy  equa- 
tion in  a model  which  postulates  a mean  flow  insta- 
bility as  the  mechanism  for  deepening.  Most  recent 
efforts,  however,  have  involved  modeling  of  the  terms 
of  the  turbulent  part  of  the  kinetic  energy  budget.  It  is 
recognized  here  that  the  Pollard  ft  al.  model  does  predict 
a possible  cyclical  steady  state.  However,  their  model 
fails  to  consider  the  turbulence  generated  above  the 
entrainment  zone  as  a source  of  energy  for  mixing 
■within  the  zone.  A mean  flow  instability  does  not  seem 
to  be  the  dominant  mechanism  for  significant  layer 
deepening,  and  simulations  using  this  model  fall  short 
of  predicting  the  observed  amount  of  deepening. 

2)  The  stable  regime  (//*> 0)  for  the  turbulent 
boundary  layer  is  not  well  understood,  especially  in  the 
limiting  case  of  retreat,  i.e.,  dh  dt^O.  Retreat  occurs 
when  the  vertical  component  of  turbulence  is  insufficient 
to  transport  heat,  momentum  and  turbulence  to  an 
earlier-established  depth  of  mixing.  Knowledge  of  the 


distribution  of  turbulent  energy  between  horizontal  and 
vertical  components  is  therefore  crucial  in  predicting 
cessation  of  mixing,  i.e.,  layer  retreat.  The  earlier  users 
of  an  equation  for  the  total  turbulent  energy  have 
neglected  this  distributional  factor. 

3)  In  most  previous  models,  all  buoyant  production 
of  turbulent  energy  has  been  consigned  to  potential 
energy  increase,  or  r=  1,  where  r is  the  ratio  of  down- 
ward entrainment  buoy-ancv  flux  to  upward  surface 
buoyancy  flux,  — bw(—h)/bu'(0).  This  could  only  lie 
possible  if  none  of  the  convectively-produced  turbulent 
energy  were  dissipated.  Numerous2  measurements  show 
that  r is  much  less  than  unity.  This  problem  is  not 
solved  by  taking  dissipation  to  be  a fixed  fraction  of 
shear  production  plus  buoyant  production  (less  buoyant 
damping)  because  layer  retreat  will  then  lie  predicted 
only  if  buoyant  damping  equals  shear  production, 
making  the  value  assigned  to  dissipation  be  equal  to 
zero.  However,  if  there  is  turbulence  available  for 
buoyancy  flux,  there  must  also  lie  dissipation.  Both 
this  apparent  dilemma  and  the  retreat  problem  are 
part  of  a single  larger  problem : dh  dt  should  not  be 
the  direct  consequence  of  an  integral  constraint  on  the 
total  turbulent  energy  equation. 

In  an  attempt  to  resolve  the  above  three  problems, 
this  paper  offers  a new  one-dimensional  model  of  the 
mixed  layer  in  which  several  processes  are  parameter- 
ized more  explicitly  than  in  previous  models.  First,  the 
fraction  of  wind-generated  turbulent  kinetic  energy- 
available  for  mixing  is  dependent  on  the  ratio  of  the 
mixed  layer  depth  to  the  Obukhov  mixing  length. 

2 Stull  (1975)  lists  experimental  observations  of  r (his  A i),  and 
the  median  value  is  between  0.1  and  0.3. 
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Second,  viscous  dissipation  is  dependent  on  a local 
Rossby  number.  Third,  separate  vertical  and  horizontal 
equations  for  turbulent  kinetic  energy  are  used,  allowing 
for  a more  explicit  treatment  of  the  mixing  process. 

Mean  turbulent  field  modeling  of  the  terms  of  these 
component  equations  is  necessary  to  put  the  model  in 
closed  form.  There  is  no  particular  precedent  for  doing 
this  in  a bulk  model,  but  Bradshaw  (1972),  Mellor  and 
Herring  (1973)  and  Lumley  and  Khajeh-Nouri  (1974) 
provide  a general  background  for  the  technique. 

In  the  model  to  be  presented  here,  conservation  of 
buoyancy  is  employed  as  a generalization  of  the 
conservation  of  heat  alone.  The  buoyancy  equation  is 
generated  from  the  heat  and  salt  equations  together 
with  an  equation  of  state, 

p=po[l—a(0—0o)+0(S—So)],  (1) 

and  the  definition  for  buoyancy, 

*=g(po— p)/po-  (2) 

In  (1)  and  (2)  6 is  temperature,  S salinity  and  p density, 
while  a and  P are  the  expansion  coefficients  for  heat  and 


salt,  respectively,  and  g is  gravity.  The  tilde  represents 
the  total  instantaneous  value  and  the  subscript  zero 
denotes  a representative  but  arbitrary  value.  The 
generalization  of  using  b rather  than  6 will  cast  the 
model  equations  in  a form  equally  applicable  to  those 
situations  where  evaporation  and  precipitation  con- 
tribute significantly  to  the  surface  buoyancy  flux  and 
the  structure  of  the  evolving  pycnocline.  The  buoyancy 
equation  also  has  a more  obvious  and  direct  role  in  the 
mechanical  energy  budget,  as  shown  in  Fig.  2. 

2.  Entrainment 

a.  Conceptual  model 

A physically  plausible  and  accurate  model  for 
predicting  the  rate  of  deepening  (or  retreat)  is  de- 
pendent upon  an  understanding  of  the  dynamics  of  the 
entrainment  process.  The  assumption  is  that  the  tur- 
bulence of  the  overlying  mixed  layer  provides  the 
energy  needed  to  destabilize  and  erode  the  underlying 
stable  water  mass.  Therefore  the  turbulent  kinetic 
energy  budget  is  the  basis  for  the  entrainment 


1 • U(0) 
work  hv 
wind  stress 

~] 

j MEAN  K.K. 

I '-A’n  (u2+v2)dz 

1 -li-l 
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plus  absorption  of  solar  radiation 


interior  buoyancy  flux 
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Fio.  2.  Mechanical  energy  budget  for  the  ocean  mixed  layer.  Asterisks  indicate  those 
processes  that  must  be  parameterized  to  close  the  system  of  equations. 
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1 d(«,+s,+i»J) 

2 dt 


-—dU  —dV-i  _ 

uw |-i rw — l+bui 

dz  dz  J 


df  /rf+tP+u?  p 


1T 


,0, 


(3) 


where  t is  viscous  dissipation  and  ensemble  mean  and 
fluctuating  components  are  denoted  by  the  upper  and 
lower  cases  respectively.  For  example,  fi=t/(z,<) 
where  U=u.  If  the  terms  of  (3)  are  known 
throughout  the  boundary  layer  then  the  evolution  of 
the  potential  energy  and  density  profile  may  be  evalu- 
ated by  using  the  budget  for  mechanical  energy. 

In  the  conceptual  model,  the  entrainment  zone 
(region  II  of  Fig.  t)  is  a region  which  is  intermittently 
turbulent  in  comparison  with  the  overlying  region  I. 
Local  turbulent  intensity  and  the  work  rate  by  this 
turbulence  on  the  interface  (z=  —h)  in  the  form  of 
—bw(—h)  is  dependent  upon  the  rate  of  supply  of 
energy  from  above,  — (d/  dz)[_u>(E/2+  p/ po)~]-k,  where 
E=  u2+t?fi-  up.  Without  this  extra  energy  source,  the 
region  will  remain  dynamically  stable,  with  a flux 
Richardson  number 


Ri=bw/(uwdU/dz-{-vwdV/dz)  (4) 


most  like!/  to  be  the  dominant  mechanism  leading  to 
observed  rates  of  entrainment. 

The  specific  mechanism  that  is  envisioned  in  the 
destabilization  of  the  interface  and  the  resulting  en- 
trainment is  a “local”  K-H  instability.  The  onset  of 
this  instability  and  its  exponential  growth  rate  is 
predicted  by  linear  two-dimensional  wave  theory.  As 
individual  wave  packets  achieve  a significant  amplitude, 
the  nonlinear  three-dimensional  effects  of  the  turbulence 
field  are  assumed  to  prevail  and  to  advect  parts  of  the 
exposed  cusps  of  denser  water  up  into  the  mixed  layer. 
Therefore,  it  is  only  the  initial  stages  of  the  instability 
that  are  strictly  of  the  K-H  type,  where  the  induced 
suction  at  the  crests  of  a perturbation  wave  on  the 
interface  is  large  enough  to  overcome  the  restoring 
buoyancy  force.  The  shear  needed  to  trigger  such  an 
instability  is  provided  by  the  local  turbulent  eddies. 
The  mean  shear  contributes  to  the  instability  but 
usually  cannot  in  itself  generate  a critical  Richardson 
number  or  directly  influence  the  frequency  and  magni- 
tude of  the  destabilizing  events. 

c.  Entrainment  hypothesis  derived 

Since  the  presumption  is  that  the  convergence  of  flux 
of  turbulent  energy  at  the  interface  is  primarily  respon- 
sible for  the  entrainment  buoyancy  flux,  the  problem  is 
to  estimate  the  time  scale  re  required  to  transport  some 
of  the  available  turbulent  energy  ( E ) to  'he  vicinity  of 
the  entraining  interface : 


larger  than  the  critical  value  for  a return  to  laminar 
flow. 

b.  The  specific  mechanism 

Of  Benjamin’s  (1963)  three  basic  types  of  instabilities, 
two  may  be  possible  here.  At  the  interface  between  the 
mixed  layer  and  the  denser  water  beneath,  a so-called 
class  A instability  will  arise  if 


where  k is  the  wavenumber  of  the  interfacial  disturbance 
and  Am  the  total  velocity  change  across  the  interface. 
From  Lamb  (1932),  the  Kelvin-Helmholtz  (K-H)  in- 
stability (Benjamin’s  class  C)  requires  a larger  value 
for  Am  : 


where  the  angle  braces  denote  the  vertical  mean  over 
the  mixed  layer,  i.e., 


(£(0)=— — f EM  dz.  (8) 

h-\-6J  -h-s 


The  mixed  layer  depth  h or  a length  scale  proportional 
to  h is  the  distance  over  which  turbulent  energy  must 
be  transported  by  the  vertical  component  of  turbulent 
velocity  tv.  Therefore  t,  is  taken  to  be  proportional 
to  h divided  by  the  rms  vertical  velocity  scale,  (it)2)1, 
giving 

T.=aiA(w*)-*.  (9) 

The  width  4 of  region  II  is  assumed  to  adjust  so  as 
to  maintain  stability : 


However,  the  class  A instability  is  dependent  upon 
energy  dissipation  in  the  lower  fluid,  and  this  is  likely 
to  be  small  compared  with  inferred  rates  of  convergence 
of  energy  flux  at  the  interface  (9/  dz)[v>(E/2-\-  p/ p0)']-h- 
For  geophysical  flows  of  this  type  having  large  Reynolds 
and  P 6c let  numbers,  the  class  C instability  is  therefore 


SAB 

RS  = Ri(-h-S<z<  -A)« 

(Al/)*+(AF)* 

= constant,  (10) 

where  A denotes  the  drop  in  the  mean  variable  across 


459 


May  1977 


ROLAND  W.  GARWOOD,  JR. 


DEPTH  AVERAGED  RICHARDSON  NUMBER.  15.8-26-2  METERS 


Fig.  3.  From  Halpern  (1974)  showing  a storm  on  20  August  that  deepened  the  mixed  layer  to  about  24  m. 


the  interface.  A too-sharp  interface  (i~0)  would  be 
dynamically  unstable  for  all  wavenumbers  ( k<S~l ) 
with  any  finite  mean  velocity  drop  | At/,  | . The  resultant 
vertical  mixing  would  increase  & until  (10)  was  satisfied. 
On  the  other  hand,  interfacial  instabilities  and  sub- 
sequent mixing  will  sharpen  (3— >0)  the  interface, 
explaining  the  very  sharp  interfaces  observed  in  grid- 
stirred  experiments  where  At/<~0.  The  combined  effect 
of  these  two  tendencies  is  to  maintain  an  equilibrium 
value  for  5 so  that  R S remains  constant.  It  is  important 
to  realize  here  that  the  assumption  that  RS  is  constant 
does  not  constitute  closure  because  the  value  of  i is  not 
a known  quantity.  This  concept  of  the  interface 
dynamics  and  the  role  of  Ri  follows  closely  the  argu- 
ment of  Csanady  (1974).  The  closure  hypothesis  of 
Pollard  el  al.,  (1973),  hA H; (AI/*+  A F2)  = constant,  is 
derived  from  (10)  only  by  making  the  additional 
assumption  that  h/ 5=  constant.  However,  the  entrain- 
ment with  shear  by  Moore  and  Long  (1971)  supports 
(10)  and  indicates  that  h 5 is  not  constant.  Halpem’s 
(1974)  measurements  of  depth-averaged  gradient 
Richardson  number  versus  time  (Fig.  3)  also  lend 
support  to  (10).  As  the  vertical  position  of  the  interface 
was  modulated  by  tidal-frequency  internal  waves,  the 
position  of  the  current  meters  may  have  passed  into  or 
through  the  interface  region,  depending  upon  both  the 
average  mixed-layer  depth  and  the  amplitude  of  the 


internal  wave.  The  storm  on  20  August  deepened  the 
layer  sufficiently  to  influence  the  envelope  for  15.8— 
26.2  m but  not  for  26.2-46.2  m. 

With  Fig.  1 in  mind,  Rf  may  be  expressed  in  terms 
of  the  bulk  model  properties  by  integrating  the  mean 
buoyancy  and  momentum  equations 


dB 

dbiv  ag 

dl 

dz  poC, 

dU 

duw 
= /F-— , 
dz 

dl 

dV 

dl 

dvw 

=-fu , 

dz 

(ID 

(12a) 

(12b) 


across  the  entrainment  zone,  from  z=  -h—i  to  z=  —h. 
If  negligible  amounts  of  momentum  and  buoyancy  are 
transported  below  the  entrainment  zone,  and  the 
interface  doesn’t  change  significantly,  then  this  inte- 
gration yields  the  so-called  jump  conditions  for  tur- 
bulent fluxes  at  the  bottom  of  the  mixed  layer: 


_ 9k 

-M-*)  = A/J-  , (13) 

dl 





- - 
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_ dh 

—uw(—h)  = i\U — , (14a) 

dt 

_ dh 

— iw(— h)=AV — . (14b) 

dt 

Using  (10)  and  the  above  jump  conditions,  one  finds 
that  shear  production  is  a fixed  fraction  of  buoyant 
damping  in  the  entrainment  zone : 

dU  • 

~^u^-(-h)=-(Kh)-'bw(-h),  (IS) 

dz 

where  (Ri)~*  becomes  the  constant  of  proportionality. 
Notice  that  within  the  accuracy  of  (13)  and  (14),  the 
flux  and  gradient  Richardson  numbers  are  equivalent 
(and  equal  to  RS)  in  the  entrainment  zone. 

Tennekes  (1973)  assumed  that  dissipation  is  an  in- 
significant part  of  the  turbulent  kinetic  energy  budget 
in  the  entrainment  zone,  implying  that  RS  must  have  a 
value  larger  than  unity.  This  zone  may  indeed  maintain 
a flux  Richardson  number  larger  than  1 and  still  sustain 
active  entrainment  because  shear  production  is  not  an 
important  local  source  of  energy  for  mixing.  Within  an 
active  entrainment  zone,  the  most  significant  source  is 
the  convergence  of  flux  of  turbulent  energy,  — (d/dz) 
Y,[_w(E/2+p/p0)~\.  Therefore  the  critical  parameter 
determining  the  rate  of  entrainment  is  not  Rf  but  is 
instead  the  ratio  P of  buoyancy  flux  to  convergence  of 
energy  flux,  i.e., 

<16) 

Then  (7),  (9)  and  (16)  give  an  entrainment  equation, 

P(-k)  = hfou(-h)/&)'(£)=tnt,  (17) 

if  P(—h)  is  the  critical  constant  parameter,  assigned 
the  value  m4. 

If  dissipation  in  the  entrainment  zone  is  either 
negligible  or  is  a fixed  fraction  of  the  flux  convergence, 

i.e., 

*(-*)-«*<£)/ T„  (18) 

then  (17)  is  also  the  consequence  of  (3),  (7),  (9),  (15) 
and  (18)  where  m,  reflects  the  combination  of  the 
constants  of  proportionality: 

(l-o,  )RS 

m,= . 

o,(R«-l) 

Eq.  (17)  is  similar  to  the  form  in  Tennekes  (1973),  with 
the  primary  difference  being  the  use  of  (u^)i(E)  rather 
than  simply  ( E )*.  This  feature  should  generalize  the 
applicability  of  the  final  mixed  layer  model,  enabling 
its  use  under  a wide  span  of  conditions,  i.e.,  the  diurnal 


and  annual  ranges  of  mixed  layer  stability  for  which 
the  ratio  (u?)/(E)  would  be  expected  to  vary  signifi- 
cantly. This  approach  to  the  entrainment  problem  is 
very  different  theoretically  from  that  of  Mellor  and 
Durbin  (1975).  Although  both  theories  employ  mean- 
turbulent-field  modeling  techniques  in  the  turbulent 
energy  budget,  Mellor  and  Durbin  neglect  the  flux 
convergence  term  altogether  and  rely  upon  a critical 
gradient  Richardson  number  at  the  interface  coupled 
with  gradient  diffusion  below  the  interface.  The  deriva- 
tion of  (17),  from  Garwood  (1976),  is  based  upon  the 
assumption  that  the  instability  mechanism  is  dynamical 
in  nature,  but  the  Mellor  and  Durbin  parameterization 
depends  upon  the  viscous  (ej^O)  character  of  the  inter- 
face, as  in  Benjamin’s  (1963)  class  A type  instability. 
Measurements  of  the  flux  convergence  term  at  an 
entraining  interface  are  nonexistent  but  are  needed  to 
settle  this  very  fundamental  difference. 

Eq.  (17)  does  not  close  the  problem  of  predicting  the 
evolution  of  the  density  structure  of  the  upper  ocean, 
given  initial  conditions,  and  surface  boundary  condi- 
tions (wind  stress  and  buoyancy  flux),  because  two  new 
unknowns,  (E)  and  (te2),  have  been  introduced. 

3.  The  bulk  equations 

Final  closure  is  achieved  with  mean-turbulent-field 
modeling  of  the  vertically  integrated  equations  for  the 
individual  turbulent  kinetic  energy  components,  plus 
the  inclusion  of  the  bulk  buoyancy  and  momentum 
equations. 

With  rectangular  coordinate  axes  having  x positive 
to  the  east,  y positive  to  the  north  and  z positive  up- 
ward, the  individual  turbulent  kinetic  energy  budgets 
are  given  by 


1 du2  dU  d /i vu2\  p du  

= — UIV — ( H 1" fljto 

2 dt  dz  dz\  2 / po  dx 

—il2uu/—t/3,  (19a) 

1 dv1  dV  d fwv2\  p dv  

= — vw 1 — H flam'  — e/3,  (19b) 

2 dt  dz  dz\  2 / po  dy 


1 dw2  d fiv*  wp\  p dw  

— bnv ( — | )-| (-  Q2uw — e/3 . (19c) 

2 dt  dz\  2 po  / po  dz 

In  deriving  these  equations,  horizontal  homogeneity 
was  assumed,  neglecting  those  terms  containing  mean 
horizontal  gradients.  The  sum  of  these  three  component 
equations  gives  (3),  the  total  turbulent  kinetic  energy 
budget.  Each  of  the  component  equations  is  usually  in 
a quasi-steady  state  because  the  dissipation  time  scale 
is  usually  much  smaller  than  the  time  scale  of  the 
surface  fluxes.  A steady-state  assumption  will  simplify 
final  solution,  but  it  is  not  actually  necessary  to  achieve 
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closure.  Notice  that  the  redistribution  terms  associated 
with  rotation,  ii„nua,  and  pressure  interaction, 


p du„ 
Po  dxa 


vanish  in  the  summation.  These  terms  therefore  did  not 
appear  and  therefore  played  no  role  in  the  earlier 
models  in  which  the  turbulent  energy  budget  was  not 
separated  into  components.  Vertical  integration  of  (11) 
and  (12a,b)  yields  the  bulk  relationships  for  mean 
buoyancy 


d{b) 

h 

dt 


dh  agQo  — 
+Ab — A bw(0), 

dt  poCp 


and  mean  momentum 


(20) 


d(U)  dh  _ 

h \-AU— a = /A<K)-«w(0),  (21a) 

dt  dt 

d(V)  dh  _ 

k +AF— A=-fh(U)-vw(0).  (21b) 

dt  dt 


Three  assumptions  were  employed  in  this  integration: 

(i)  Vertical  fluxes  are  negligible  below  the  mixed 
layer.  Therefore, 


lrw(—h—S)  = uw(—h  — i)  = vw(—h—6)=0. 

(ii)  The  mixed  layer  is  sufficiently  homogeneous  so 
that 

A /**<«)-fl(-A-8), 

A U*>(U)-U(-h-S), 

A V~(V)-V(-h-i). 

(iii)  Horizontal  homogeneity  is  assumed  for  all  mean 
variables.  The  approximation  of  local  horizontal 
homogeneity  with  regard  to  the  mean  turbulence  fields 
is  usually  an  accurate  assumption  because  of  the  short 
time  scale  for  the  turbulence.  On  the  other  hand,  it  is 
recognized  that  the  mean  buoyancy  and  momentum 
fields  are  not  one-dimensional  for  all  time  and  space 
scales.  However,  such  advective  effects  shall  be 
neglected  here  in  order  to  emphasize  some  new  aspects 
which  are  fundamental  to  the  one-dimensional  model. 


4.  Closure  hypotheses 

Mean  turbulent  field  modeling  of  the  terms  of 
(19a-c),  integrated  across  the  mixed  layer,  will  close 
the  problem.  That  is,  there  will  be  five  equations  for 
the  five  unknowns 

h,  ( B >,  <0-<l/>+i<V>,  (w*>,  (i5+i»). 

In  addition,  there  will  be  empirical  constants  to  be 
specified. 


a.  Viscous  dissipation 

For  fully  turbulent  geophysical  flows  having  large 
Reynolds  numbers,  viscous  dissipation  of  the  turbulence 
occurs  primarily  in  the  small  eddies  which  are  locally 
isotropic.  As  explained  by  Tennekes  and  Lumley  (1972), 
an  estimate  of  dissipation  is  made  by  taking  the  rate  at 
which  the  largest  eddies  supply  energy  to  the  smaller 
eddies  (equal  to  the  rate  of  dissipation)  to  be  propor- 
tional to  the  reciprocal  of  the  time  scale  of  the  largest 
eddies.  The  net  rate  of  dissipation,3 


r°  du,  du, 

D=  tdz*=  „ dz,  (22) 

J -K-i  J -K-i  dxj  dxj 

and  the  vertical  mean  of  turbulent  energy,  (£),  are 
accordingly  used  to  define  a dissipation  time  scale, 
t. = <£)/(«),  or 


(23) 


1)  Dissipation  in  shallow  mixed  layers,  Ro»1 

If  the  time  scale  of  these  largest  eddies  is  proportional 
to  the  mixed  layer  depth  divided  by  the  rms  turbulent 
velocity,  i.e., 

t,= /><£)-»,  (24) 

then  an  integral  model  for  dissipation  in  the  mixed 
layer,  independent  of  viscosity  and  the  small  scales,  is 


tdz  = mi(E),t 


(25) 


where  f»i  is  a constant  of  proportionality.  For  those 
situations  where  the  turbulent  velocity  scale  (£)*  is 
proportional  to  the  surface  friction  velocity  Eq.  (25) 
is  equivalent  to  the  parameterization  used  by 
Miropoi’skiy  (1970)  and  Denman  (1973).  Such  is  the 
case  only  for  neutral,  fw(0)  = 0,  mixed  layers. 


2)  A LIMITING  DISSIPATION  TIME  SCALE  FOR  DEEPER 
(Ro~l)  MIXED  LAYERS 

In  deeper  boundary  layers,  planetary  rotation  turns 
the  mean  shear  direction  with  depth  and  thus  influences 
the  geometrical  aspects  of  the  integral  scale.  This 
introduces  a second  integral  time  scale 

W"1,  (26) 

the  inverse  of  the  coriolis  parameter. 

It  is  becoming  increasingly  clear  from  such  studies 
as  Arya  and  Wyngaard  (1975)  and  Sundararajan  (1975) 
that  this  rotational  time  scale  plays  an  important  role 
in  the  internal  structure  of  the  convective  planetary 

* This  D is  equivalent  to  the  7>*/po  of  Kraus  and  Turner  (1967), 
where  D*  is  the  total  rate  of  dissipation  per  unit  area. 
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boundary  layer  or  mixed  layer.  The  concern  here  is  more 
with  the  bulk  properties  of  the  region  and  less  with  the 
details  within  the  mixed  layer,  but  it  is  suggested  that 
this  time  scale  has  an  important  role  in  the  overall 
turbulent  energy  budget.  The  mean  shear  profile  and 
the  turbulent  energy  budget  are  inseparable  because 
of  the  link  through  local  shear  production. 

Rather  than  to  simply  replace  the  convective  scale 
ti  by  the  rotational  scale  r2,  it  is  assumed  that 

(r.)-1=(r1)-‘+(r1)-‘.  (27) 

This  is  the  simplest  combination  of  the  two  scales  that 
retains  the  asymptotic  characteristic  of  r,  —*  tj  as 
hf  — * 0.  This  gives 


c.  Shear  production 

The  vertical  integral  of  shear  production  reduces  to 


-fj 


—dU  _dF\ 


zuJSU  (0), 


where  SU (z)  is  the  magnitude  of  the  mean  velocity 
associated  with  mean  shear,  i.e., 

&U(z)  = tlP+V1-(lP)-(Vt)y.  (33) 

In  this  instance,  the  inhomogeneity  of  the  mean  velocity 
field  cannot  be  neglected.  An  additional  source/sink 
term  is 


D=  [ «/z=m,(£)»+«,/A(£)  (28)  dJ_W(Po  + 2 ) KPo+2  )]0 

J -k-l 


£>=mi(£)*^l 


m u*  \ 

m,  (£)V 


PO  J-A 


«* 

Ro= — (29) 

*/ 

is  a Rossby  number  for  the  turbulent  boundary  layer. 

b.  Redistribution  of  turbulent  energy 

The  vertical  integral  of  the  pressure  redistribution 
term, 

r°  p dua 

Ra=  I —dz,  (30) 

J -*-»  PO  bxa 

is  an  important  source  or  sink  term  for  the  individual 
turbulent  kinetic  energy  budgets,  even  though 
•Ri+£2+£i— 0. 

Following  the  early  lead  of  Rotta  (1951),  and  in 
agreement  with  the  dominant  term  of  the  rational 
closure  technique  of  Lumley  and  Khajeh-Nouri  (1974), 
the  bulk  formulation  is 

£«=m2(£>»«£)-3<^)).  (31) 

In  addition  to  dimensional  consistency,  the  concept 
leading  to  (31)  is  that  of  a “return  to  isotropy.”  In 
other  words,  the  pressure-strain  rate  interaction  tends 
to  restore  equal  distribution  of  energy  among  the  three 
components.  This  interaction  is  expected  to  be  some- 
what dependent  upon  stability,  but  this  is  assumed  to 
be  a higher  order  effect  and  is  neglected  here.  The 
rotational  redistribution  terms  il.uu„  are  also  assumed 
to  be  of  a higher  order  and  are  neglected  at  this  point. 
Their  inclusion  would,  however,  create  the  intriguing 
possibility  of  entrainment  rate  being  susceptible  to 
wind  direction. 


The  surface  term  reflects  a source  attributable  to  break- 
ing surface  waves,  and  it  might  be  modeled  as  being 
proportional  to  «,*  for  a fully  developed  wave  field.  The 
term  at  z=  —k—S  accounts  for  a possible  loss  due  to 
downward-radiating  internal  wave  energy.  Stull  (1975) 
found  that  this  may  be  significant  in  a rapidly  deepening 
unstable  atmospheric  boundary  layer.  The  term  will  be 
neglected  here. 

If  41/(0)  is  proportional  to  «*,  then  (32)  may  be 
combined  with  (34)  to  give  net  “wind-generated”  rate 
of  production4 : 

T — dU  dV  d /wp  w£\-i 

G — I — uw 1 -vw 1 — ( 1 ] viz 

J -K-i  L dz  dz  dz\  po  2 / J 

(Auy+{Avy  ah 

=m2w*,-| . (35) 

2 dt 

The  last  term  of  (35)  comes  from  integrating  the  shear 
production  across  the  entraining  interface  (Niiler, 
1975). 

d.  Summary  of  modeled  equations 
A final  set  of  equations  has  been  generated : 

The  extrainment  buoyancy  flux,  from  (17)  is 

— . ,,  m^vpyiE)  (36) 

—bw(—h)  = . 

h 

The  budget  for  the  horizontal  components  of  turbulent 
kinetic  energy  comes  from  (35a),  (31),  (28)  and  the 


4 This  G is  equivalent  to  Kraus  and  Turner’s  G*/po  plus  the 
term  added  by  Niiler. 
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sum  of  (19a)  and  (19b)  vertically  integrated,  giving 


5.  General  behavior  of  the  equations 


(hiu'+v1))  =*»,«,*- 

2 at 


bw{— A)  | AC7 1 * 
2AB 


-w,«£)-3<w*))<£)»- 


2mi/ 


, (£)»+ 


m6  \ 

-/*)<£>• 

m i / 


(37) 


The  budget  for  the  vertical  component  of  turbulent 
kinetic  energy  comes  from  (31),  (28)  and  the  vertical 
integral  of  (19c),  giving 

1 d _ _ 

(A(w*))  = \hbw  ( — h)  — \hutbt  +m2((E)~i  (w1))  (£)* 

2 a; 


<£)»+-/* 


mi 


(38) 


The  mean  buoyancy  and  complex  mean  momentum 
equations  are 


d(B) 

— , — as 

t 

—bw{ — h) — bu>(0)+ Qo, 

(20) 

at 

poCp 

d(C) 

P 

at 

= cw(—h)—cw(0)—if{C)k. 

(21) 

The  jump  conditions  relate  the  entrainment  fluxes  to 
the  mean  momentum  and  mean  buoyancy  and  the 
rate  of  deepening, 


_ dh 

-cw(-A)=AC— , (14) 

at 

_ dh 

— 6w(— A)  = AB — . (13) 

at 

The  values  of  uw(0),  iro(O),  bw(0)  and  the  radiation 
absorption  Q(z)  are  assumed  to  be  given  time-dependent 
variables  from  which  are  determined 

«**=  | ttif(0)-Hiw(0)  | , (39) 

u»bt—~ bw(0) —[  (q f Qdx)dz.  (40) 

poCV7_*_j  \ h-\-Sj , ) 

Also  assumed  to  be  given  are  the  mean  buoyancy  and 
momentum  just  below  the  mixed  layer:  B(—h—b) 
and  C(—  h— b).  Therefore 

AB  = (B)-B{~h-b)  l 
AC=(C)— C(-A-3)  I' 


a.  N bndimensional  form  of  the  turbulent  energy  and 
entrainment  equations 

Using  the  given  flux  scales  u,  and  6*,  new  dimension- 
less variables  are  defined : 


H*  = 


utb*h 

2 WjM*3 


hAB  dh 
P*= 


2wj«»3  at 

/miy  (E) 
\wj  M*2 

©'  (u>2) 

u* 


—hbw{—h) 

2msuJ 


(41) 

(42) 

(43) 

(44) 


The  mixed  layer  stability  parameter  H*  is  the  ratio  of 
the  effective  buoyancy  flux  (from  surface  heating  and 
net  precipitation  minus  evaporation)  to  wind-stress 
production.  It  is  proportional  to  h/ L where  L is  the 
Obukhov  length.  Both  L and  H*  may  be  negative 
should  there  be  a positive  surface  buoyancy  flux  due 
to  surface  cooling  and/or  a net  evaporation  minus 
precipitation.  A value  of  zero  for  H*  (L— >«)  repre- 
sents a neutral  mixed  layer.  The  interface  stability 
parameter  P*  is  a dimensionless  measure  of  the  entrain- 
ment rate.  The  parameters  Efti  and  Ef33  are  measures  of 
the  total  turbulent  kinetic  energy  and  the  vertical 
component  of  turbulent  kinetic  energy  respectively. 
The  values  of  P*  and  Efaa  will  depend  upon  the  values 
of  H*  and  Ro-1  that  are  determined  by  the  prescribed 
values  of  «»2  and  together  with  the  current  value 
of  h. 

Invoking  the  quasi-steady-state  assumption  for  the 
turbulent  energy  budget,  the  entrainment  and  turbulent 
energy  equations  (36),  (37)  and  (38),  become 

P*=iplE!,(E‘3 ,)»,  (45) 

P* 

0 = 1 +— pt  (£*  -3£j3)  (£*  ) * 


-|£M+/»Ro-‘),  (46) 
0=  ~H*~P*+p3{Efu-iBr33)(Pru)\ 

-K,  (£?/+/>,  Ro-1),  (47) 

where 


Ri*=AAB/|AC|2, 


pi=mt/mi, 

(48a) 

Pt= 

(48b) 

pt=mi/mi. 

(48c) 
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from  the  Kato  and  Phillips  (1969)  laboratory  results, 
dh 

hAB—  (6*=0)  = 2.5w*3,  (54) 

dt 

together  with  P*(H*=  0)  from  Fig.  4, 
hAB  dh 

P*(H*  = 0)  = (6,  = 0)  = 0. 165, 

2m3ut3  dt 

giving 

m,=  1.25/P*(fl*=0)=7.6.  (55) 

A check  on  this  value  comes  from  (32).  Eqs.  (32)  and 
(35),  with  Ri*>2,  give 

ni3=&U  (0  )/m*.  (56) 

Therefore,  a value  of  w3=7.6  is  quite  reasonable. 

c.  i.eneral  solution — including  dissipation  enhancement 

The  general  solution  to  P*  is  largely  a function  of 
only  two  variables,  H*  and  Ro_!  if  Ri*7>2.  This  solu- 
tion, from  Garwood  (1976),  is  given  in  Fig.  6.  It  is  the 
result  of  the  simultaneous  solutions  of  (45),  (46)  and 
(47).  The  earlier  solution  (Fig.  4)  is  represented  by  the 
limiting  case,  P*(Ro~'=0).  In  examining  Fig.  6, 

Z*=/>3 Ro-1  may  be  regarded  as  the  nondimensional 
mixed-layer  depth  and  H*  may  be  considered  a measure 
of  stability.  Two  new  features  appear  in  the  general 
solution.  First,  the  rate  of  entrainment  decreases  with 
increased  Z*.  Second,  in  the  retreat  mode  (P*=0), 
H’mtt  is  not  a constant  but  is  a function  of  layer  depth, 
rotation  and  wind  stress.  In  particular,  a neutral  steady 
state,  P*—0  for  //*= 0,  is  predicted  for  a sufficiently 
large  value  of  Z*.  However,  starting  with  A««*/ /, 
neutral  steady  state  would  be  approached  only  after  a 
relatively  long  time.  This  is  not  likely  to  be  achieved  in 


TIME  (DAYS) 

2 

diurnal  heating 
cycle 


no  surface 
heating 


Fig.  7.  Mixed  layer  depth  versus  time  for  two  h\  pothetical  cases 
having  constant  wind  stress  and  an  initial  linear  siratin - a'.  ■ Tin- 
case  with  a diurnal  heating  cycle  exhibit-  a modulation  : the 
longer  term  rate  of  deepening,  as  predicted  In  ~ 


geophysical  flows,  except  perhaps  by  the  .Mmi-nherii 
boundary  layer  during  the  polar  wintet  !’■  at.  : 

Arya,  1974). 

6.  Specific  differences  in  comparison  to  earlier 
models 

a.  Nonlinearity  in  P*(H*)  in  deepening  mixed  liver.' 

This  nonlinearity  is  present  with  or  without  rotation. 
Therefore  it  is  sufficient  to  examine  the  vi  .pier  case 
with  Ro7i>l  (Fig.  4)  to  demonstrate  this  effect.  If  II*  is 
perturbed  by  a fluctuating  surface  buoyant  y flux  such 
as  that  associated  with  the  diurnal-period  heating 
cycle,  then  the  mean  P*  over  a complete  cycle  will  be 
le4s  than  P*  for  the  mean  //*,  i.e., 


P*{H*)^P*(H*).  (57) 

Fig.  7 demonstrates  the  consequence  of  this  phe- 
nomenon over  the  course  of  a few  days.  This  non- 
linearity will  result  in  the  modulation  of  the  long-term 
trend  by  shorter  period  fluctuations  in  stability, 
particularly  those  of  the  diurnal  cycle. 


. i§§§ m H X 


V ; : V-V\  V->'  \'T  . >■ 

n Mm  i ! 

Fig.  6.  General  solution  to  the  entrainment  and  turbulent 
kinetic  energy  equations  (45),  (46)  and  (47),  if  Ri*»l  and 


b.  Cyclical  steady  state 

For  the  sake  of  studying  the  relative  response  of  the 
mixed  layer,  all  cycles  will  be  assumed  to  be  repetitive 
(in  steady  state).  A sinusoidal  surface  buoyancy  flux5 
and  a constant  wind  stress  are  used  to  drive  the  model, 
i.e.,  |cw(0) | =tt*2=  constant  and  —bic(Q)= 

Xexp(— iwt).  In  other  words,  in  this  hypothetical  case, 
bw( 0)  is  assumed  to  have  no  mean  components.  The 
period  2ir/w  is  not  specified,  but  the  results  will  be 
particularly  relevant  to  the  annual  heating/cooling 
cycle.  The  diurnal  cycle  is  not  as  likely  to  be  in  a steady 
state,  but  the  result  will  also  approximate  that  response. 


‘The  effective  surface  buoyancy  flux  is,  for  the  sake  of  con- 
venience, assumed  here  to  include  the  solar  radiation  component 
(as  if  all  shortwave  radiation  were  absorbed  at  the  immediate 
surface). 
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Fig.  8.  The  cyclical  response  of  /»(/)  and  the  corresponding 
closed  loop  in  the  (2*,  H *)  plane  due  to  an  imposed  periodic 
surface  buoyancy  flux  —bw(0). 

Cyclical  steady  state  is  depicted  by  Fig.  8.  This 
shows  the  locus  of  coordinates  ( 7*,H *)  as  the  effective 
surface  buoyancy  flux  progresses  through  a heating 
[6mi(0)<0]  and  then  a cooling  [6k>(0)>0]  phase.  The 
initial  [0]  and  final  [4]  points  coincide  at  the  point  for 
the  neutral  steady  state.  On  an  annual  basis,  this  would 
also  correspond  to  the  vernal  equinox.  Between  [0] 
and  [1]  the  buoyancy  flux  is  directed  downward 
(H*> 0)  and  is  increasing,  causing  the  mixed  layer  to 
retreat  to  its  minimum  depth  at  [1].  Between  [1]]  and 
[2]  there  is  active  entrainment,  dh/dt>0,  but  the  rate 
of  deepening  is  very  slow  because  of  continued  down- 
ward surface  buoyancy  flux.  In  spite  of  this  downward 
surface  heat  flux,  the  sea  surface  temperature  will  begin 
to  drop  as  soon  as  entrainment  heat  flux  becomes 
dominant,  i.e.,  —6ui(—h)>—dw( 0).  This  will  occur 
prior  to  the  autumnal  equinox  [2].  Because  of  the  large 
seasonal  buoyancy  gradient  created  during  the  summer 
when  the  heat  added  through  the  surface  was  not 
mixed  deeply,  h does  not  increase  rapidly  until  after 
the  surface  heat  flux  maximum  at  £3].  Of  course, 
autumn  storms  (with  greatly  increased  w*3  and  en- 
hanced evaporation)  could  more  quickly  overcome  this 
buoyancy  gradient  and  accelerate  the  rate  of  deepening. 
Any  particular  oceanic  region  has  neither  a constant 
wind  stress  nor  a perfectly  sinusoidal  buoyancy  flux. 
Any  mean  (over  one  year)  buoyancy  flux  can  only  be 
accounted  for  by  including  advection,  which  has  not 
been  done  here.  However,  the  purpose  here  has  been  to 
show  that  a one-dimensional  mixed  layer  model  can 
explain  the  observed  cyclical  steady  state.  The  inclusion 


of  advection  and  more  realistic  surface  forcing  would 
alter  to  some  extent  ihe  shape  and  relative  position  of 
the  predicted  closed  locus  in  the  ( Z*,H *)  plane. 

The  relative  importance  of  the  rotational  and  surface 
buoyancy  flux  scales  needs  to  be  examined.  The  magni- 
tude of  the  mixed  layer  response  is  found  to  be  a 
function  of  the  ratio  of  the  rotational  length  scale  («,/ /) 
to  the  buoyancy  flux  scale  (w,V6*),  i-e.. 


Fig.  9 shows  the  cycle  of  h(t),  nondimensionalized  on 
u+/f  as  a function  of  B*. 

For  large  /?*,  which  is  typical  of  the  annual  cycle  in 
temperate  oceanic  regions,  the  shallow  mixed  layer 
depth  at  w(=  ir/2  is  attributable  to  the  strong  influence 
of  the  buoyant  damping  corresponding  to  a small 
positive  Obukhov  length  scale.  A more  realistic,  non- 
constant would  influence  the  time  for  the  occurrence 
of  this  minimum  layer  depth.  In  addition  to  the  change 
in  the  range  of  h with  B*,  there  is  an  effect  upon  the 
shape  of  the  function  h(t).  For  the  case  of  weak  heating 
and  cooling,  the  variation  in  h is  almost  in  phase  with 
bw( 0).  In  the  (Z*,  H*)  plane,  the  locus  of  the  cycle  is 
small  and  elongated.  For  increasingly  larger  B*,  the 
heat  and  therefore  buoyancy  is  stored  during  the 
“summer”  at  an  increasingly  shallower  depth  (in  com- 
parison to  A,hm~«*//).  This  creates  a hysteresis  effect 
by  retarding  the  subsequent  rate  ol  deepening.  This 
phase  lag  in  the  heat  storage  naturally  has  important 
implications  for  the  interaction  with  the  atmosphere  for 
all  surface  flux  time  scales,  from  one  day  to  periods  of 
climatological  importance. 


PHASE  (out) 


Fig.  9.  Cyclical  steady-state  response  to  a sinusoidal  surface 
buoyancy  flux  and  a constant  wind  stress.  Mixed  layer  depth  is 
nondimensionalized  on  «,//  and  B*=  | H.ft.  1 /(/«.*)  is  a measure 
of  the  strength  of  the  annual  heating/cooling  cycle. 
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7.  A theoretical  framework  for  model  comparison 
and  testing 

Fig.  6,  in  which  the  solution  for  entrainment  rate  is 
depicted  as  a function  of  stability  (//*)  and  a second 
parameter  Z*  (=pjRo_l  in  the  present  formulation), 
poses  a framework  for  model  intercomparison.  Fig.  10 
shows  illustrations  of  the  P*  surface  as  predicted  by 
three  representative  models. 

Fig.  10a  is  that  for  the  Kraus-Turner  (KT)  type  of 
model.  This  is  representative  of  all  models  (including 
those  of  Miropol’skiy,  Denman  and  Niiler)  for  which 
shear  production  G is  parameterized  as  being  propor- 
tional to  ut3,  and  net  dissipation  D is  a fixed  fraction 
of  G,  giving 

P^T—Ci—CiH*.  (59) 

In  this  case  there  is  no  dissipation  enhancement  with 
increased  layer  depth  or  rotation,  as  measured  by  the 
second  variable  Z*. 

Fig.  10b  is  for  the  Elsberry  et  al.  (EFT)  model  in 
which  Da  1— exp(— Z*).  As  in  the  KT  model,  P*  de- 
pends linearly  upon  H*,  but  the  rate  of  deepening  is 
checked  with  increased  Z*=  h/Z : 

Pv.rr—  Ci  exp(—Z*)—c2H*.  (60) 

This  dissipation  enhancement,  however,  is  still  in- 
sufficient to  predict  a cyclic  steady  state  because  the 
locus  of  P*=  0 fails  to  cross  the  //*=0  line. 

Fig.  10c  is  for  the  model  of  Kim  (1976)  having  a 
constant  background  dissipation  to  in  addition  to  the 
fixed  fraction  of  shear  production,  giving 

K =c1-c3Z*-c2H*.  (61) 

Notice  that  (61)  predicts  a possible  steady  state,  both 
neutral  and  cyclical,  because  the  P*—0  locus  crosses 
the  H*=  0 line  at  Z*=  ci/ci.  Kim’s  background  dissipa- 
tion is  not  necessarily  as  strong  at  that  of  EFT  for 
small  Z*,  but  this  parameterization  denotes  a stronger 
net  dissipation  with  increasingly  larger  values  of  Z*. 

If  the  dissipation  enhancement  effects  of  EFT,  Kim 
and  Garwood  are  assumed  to  be  physically  equivalent 
and  Z*=p3 Ro-1,  then  Elsberry’s  Z and  Kim’s  to  are 
both  related  to  planetary  rotation,  i.e., 

Zee  ujf,  (62) 

<oa  /«**•  (63) 

Resnyanskiy  (1975)  also  suggests  that  this  second 
length  scale  Z should  Ire  proportional  to  u,/f.  It  is  an 
interesting  development  to  find  that  present-day 
workers  are  returning  to  the  rotational  length  scale 
originally  proposed  by  Rossby  and  Montgomery  (1935). 

The  framework  P*(Z*,  //*)  is  also  applicable  to  data 
analysis.  An  empirical  determination  of  the  P*  surface 
as  a function  of  Z*  and  //*  is  needed  to  test  and  further 
improve  the  one-dimensional  model.  Actually  accom- 
plishing such  an  empirical  fit  is  difficult  because  it 
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Fig.  JO.  Solutions  to  the  entrainment  function  P*  in  the 
(//*,  Z*)  plane  for  the  models  of  (a)  Kraus  and  Turner  (1967), 
(b)  Elsberry  et  al.  (1976)  and  (c)  Kim  (1976).  Compare  these 
with  Fig.  6. 


requires  identifying  the  true  depth  of  the  turbulent 
boundary  layer.  This  may  only  occasionally  coincide 
with  the  apparent  h(l)  because  the  diurnal  variation  in 
stability  will  give  rise  to  a diurnal  retreat  and  deepening 
cycle,  not  easily  discernable  in  temperature  profiles. 
Nevertheless,  this  information  is  needed  to  adequately 
test  a model.  Any  particular  mixed  layer  model  may  be 
calibrated  to  predict  the  apparent  h(t)  at  a suitable 
site  as  long  as  the  variations  in  H*  and  Z*  are  small, 
but  such  a model  may  be  of  little  use  under  widely 
varying  conditions. 

8.  Conclusions 

This  new  model  for  the  ocean  mixed  layer  suggests 
answers  to  the  questions  raised  earlier  concerning  the 
general  applicability  of  bulk  models  based  upon  the 
turbulent  energy  budget : 

1)  Planetary  rotation  is  assumed  to  influence  the 
dissipation  time  scale  for  the  turbulence.  This  enhances 
dissipation  for  deeper  mixed  layers  and  enables  a 
cyclical  steady  state  on  an  annual  basis. 

2)  The  rate  of  entrainment  (P*)  for  the  stable 
regime,  H*>  0,  is  not  accurately  reflected  by  a linear 
extrapolation  of  the  H*^  0 situations.  This  is  particu- 
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larly  important  in  modeling  the  oceanic  boundary 
layer.  Unlike  the  atmospheric  boundary  layer  case, 
most  of  the  solar  radiation  does  not  penetrate  the  layer. 
Therefore  downward  turbulent  heat  flux  in  the  oceanic 
boundary  layer  is  as  important  as  the  upward  flux 
during  the  course  of  both  the  diurnal  and  annual  cycles. 
The  nonlinearity  of  P*(U*),  which  is  greatest  for  H*> 0 
results  in  a modulation  of  the  long-term  trend  of  h(t) 
by  the  diurnal  component  of  surface  heat  flux. 

3)  In  this  model,  buoyant  production  is  somewhat 
more  efficient  than  shear  production  as  a source  of 
energy  for  vertical  mixing  because  of  its  unique  effect 
on  the  vertical  component  of  the  turbulent  velocity. 
However,  this  efficiency  is  much  less  than  the  100% 
denoted  by  the  prototype  model. 

Finally,  a framework,  P*(Z*,  II*)  has  been  suggested 
for  model  comparison.  This  is  also  a potential  basis  for 
data  analysis  and  experiment  design.  While  the  models 
do  have  didactic  value,  knowledge  of  P*(Z*,  H*)  from 
observations  alone  would  be  sufficient  for  the  prediction 
of  layer  deepening  and  retreat. 
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